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Random walks are one of the best investigated dynamical processes on graphs. A particularly 
fascinating phenomenon is the scaling relationship of fluctuations a with the average flux (/). Here 
we analyze how network topology and nodes with finite capacity lead to deviations from a simple 
scaling law a ~ (/)^- Sources of randomness are the random walk itself (internal noise) and the 
fluctuation of the number of walkers (external noise). We obtained exact results for the extreme 
case of a star network which are indicative of the behavior of large scale systems with a broad degree 
distribution.The latter are subsequently studied using Monte Carlo simulations. We find that the 
network heterogeneity amplifies the effects of external noise. By computing the ‘effective’ scaling of 
each node we show that multiple scaling relationships can coexist in a graph with a heterogeneous 
degree distribution at an intermediate level of external noise. Finally, we analyze the effect of a finite 
capacity of nodes for random walkers and find that this also can lead to a heterogeneous scaling of 
fluctuations. 


INTRODUCTION 

Complex networks m consist of a fascinating re¬ 
search topic which has, during the last decade, revo¬ 
lutionized our understanding of dynamically interacting 
systems. Several complex computational and biological 
networks are, in fact, transport networks meaning that 
network edges serve as channels of flux towards selected 
nodes. The Internet, for example, serves as an informa¬ 
tion transport network with the edges transferring infor¬ 
mation flux, measured in bytes per second, from node 
to node. Thus, transport phenomena on complex net¬ 
works comprise a characteristic research subfield which 
has attracted many researchers as, besides its theoretical 
interest, it has also important engineering applications 

[ZHS]. 

For many complex, network-like systems the flow of 
material or information is an important functional fea¬ 
ture 0 HOI HI]. Examples include the flow of car¬ 
bon atoms (in the form of diverse chemical compounds) 
through metabolic networks nans], the flow of infor¬ 
mation through the network of internet routers mus], 
traffic flow in streets or roads [T6HT8] and via train con¬ 
nections m, the flow of material through machine net¬ 
works in industrial production [T6j [20] and many more. 

As a consequence, predicting or understanding the pat¬ 
tern of fluxes (i.e. the material/information transfer 


through all nodes in the network) in terms of the net¬ 
work topology is a principal goal in the investigation of 
dynamics on graphs. 

It has been observed, for example, that the distri¬ 
bution of metabolic fluxes in the bacterium Escherichia 
coli follows a power law, similar to the degree distribu¬ 
tion of the underlying metabolic networks [21]. Maps of 
random walks reveal the community structure in com¬ 
plex networks [22] [23]. In scale-free graphs, for exam¬ 
ple, excitations can self-organize into wave-like patterns 
around hubs m- A network derived from passenger flow 
can serve as a foundation for predicting the spread of 
epidemic diseases [25]. Clearly, many differences exist 
among the above systems e.g., on the level of conserva¬ 
tion laws, the typical signal-to-noise ratio, system size 
and the relevant architectural properties of the corre¬ 
sponding networks. 

On the theoretical level, over the last decade substan¬ 
tial progress has been made in understanding such flow 
on graphs. Elux-balance analysis, a method for predict¬ 
ing the steady-state distribution of metabolic fluxes for 
a given metabolic network from nutrient availability and 
an objective function (e.g., maximizing growth rate) is 
a highly successful tool for distinguishing between lethal 
and viable mutants in simple organisms [2i[27]. This 
method has applications as far reaching as the metabolic 
states of human cells [28]. Its’ variants have been used 
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to study which network features are enhanced during a 
simulated evolution of simple flow networks, when requir¬ 
ing robustness against link or node removal [29] or as a 
function of task complexity [30] . 

A strategy for investigating, how network topology af¬ 
fects the flux pattern, is to explore simple dynamics on 
graphs, serving as benchmarks for these investigations. 
Random walks are an important class of such benchmark 
models which are proven very useful for disentangling the 
universal network effects from those effects specific to in¬ 
dividual application systems [1ZIISIII22]. We would like 
to stress that for random walks on graphs, in addition 
to numerical simulations, also a substantial number of 
mathematical results exists (see, e.g., the review by Lo- 
vasz [33]b 

In a relatively recent article [34| several real transport 
networks were studied including the Internet, micropro¬ 
cessor networks, river streamflow networks and highway 
networks among others. The authors were interested in 
the relationship between the average number of packets 
{fi) arriving at a node i during a certain time interval i.e. 
the average node flux, and the standard deviation of 
this quantity. They find that (fi) has a broad distribu¬ 
tion and that scales as a power law with the average 
flux, i.e. 

o- ~ (/)“ (1) 

It was suggested [34| that physical systems are divided 
in roughly two classes, those with a = 1/2 and those 
with a = 1. The power law behavior with a = 1/2 is at¬ 
tributed to “internal” fluctuations while the a = 1 is re¬ 
lated to strong “external” noise. Internal noise is caused 
by the fact that the system itself consists of discrete units 
and it is inherent in the very mechanism by which the 
system evolves. External noise denotes fluctuations cre¬ 
ated in a system by the application of a random force, 
whose stochastic properties are supposed to be known. 
Subsequent studies [SDESIISS] refined this view in the 
following way m- The relation between the dispersion 
and the nodes fluxes can be separated in two parts, one 
due to “internal” and the other due to “external” noise. 

= ( 2 ) 

where the theoretically predicted value of the parameter 
c is equal to zero in the absence of external noise. 

An important question here is: When can one call data 
‘represented by / obeying a power law ’[38]? It is there¬ 
fore instructive to look in more detail at the respective 
explicatory power of the two ‘models’, the highly aggre¬ 
gated single scaling law parametrization given by Eqa 
and the overlay of two scaling relationships as represented 
by Eqj^ In any case, the main idea that the scaling re¬ 
lation between the dispersion and the nodes fluxes can 
be used as a means to study the collective dynamical 


properties of a large network and the interplay between 
“internal” and “external” noise is rather appealing and 
has triggered a substantial amount of applied research 

IMS]. 

In this paper we are particularly interested in the ef¬ 
fect of network heterogeneity on the dispersion-flux rela¬ 
tion. We use scale-free networks with degree distribution 
p{k) ^ k ^ for different 7 which allow us to model var¬ 
ious ranges of heterogeneous structures since networks 
with 7 = 2 have a large number of very highly connected 
nodes while networks with 7 = 3 have much fewer nodes 
with considerably lower maximum degree. Such networks 
are markedly different from random networks or lattices 
where all nodes have statistically the same properties and 
,hence, their topology is uniform and homogeneous. We 
use a model of multiple walkers performing fixed length 
random walks on network structures, similar to the model 
proposed in [34]. The number of walkers W may be ei¬ 
ther fixed or a random variable uniformly distributed in 
the interval [W — AW, IT + AW]. We show that the 
model is exactly solvable for the star network and this 
solution provides valuable insights for the dispersion-flux 
relation on scale-free networks with 2 < 7 < 3. Based 
on the intuition provided by the results for the star net¬ 
work we performed Monte Carlo simulations of random 
walks on scale-free networks. We confirm that the ran¬ 
dom walk model leads to dispersion-flux power law scal¬ 
ing with a = 1/2 when AW = 0, to a = 1 when AW is 
large and we are able to observe intermediate exponents 
for a capable AW spectrum. We also find, that an al¬ 
ternative analysis form (Eq. is a rather effective way 
of data analysis and we show that network heterogeneity 
enhances the transition to the a = 1 regime requiring 
lower AW for lower 7 exponents. Einally, we show that 
simple modifications of the random walk model with the 
inclusion of a node capacity or excluded volume interac¬ 
tions lead to regimes with non-power law dispersion-flux 
scaling. 

METHODS 
Random walk model 

In order to study how Eqj^ arises and its range of va¬ 
lidity we use a dynamical model initially proposed in [34] 
and subsequently used in [m [36] . We study the diffu¬ 
sion of multiple walkers starting with a network of size 
N. We use scale-free networks as the diffusion substrate 
in order to ensure a broad degree distribution. In the 
simplest version of this model, we assume that a num¬ 
ber of w walkers are randomly placed on the nodes of 
the network. This number w is the realization of a ran¬ 
dom variable chosen with uniform probability from the 
range [W — AW, W + AW]. Each node has a capacity 
C to accept walkers, where C is the number of walkers 
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FIG. 1: (A) Monte Carlo Simulation results on a 
scale-free network with 7 = 3.0, kmin = 2 and size 
N = 10000. W = 5000 walkers performed 100 
M = 1000 step walks. Lines are the ‘time series’ of the 
walker arrivals at the counters of 2 sample nodes with 
degree /c = 46, 30 respectively. (B) A star network with 
= 10 nodes. Node 0 is the central hub with degree 
ko = 9. All other nodes have degree k = 1. 


that can simultaneously be on the same node. In the 
case that C is larger than the total number of walkers 
the problem is equivalent to the diffusion of w indepen¬ 
dent walkers. When (7 = 1 the problem is equivalent 
to diffusion with excluded volume interactions. Unless 
explicitly stated otherwise the results of our simulations 
are for unlimited capacity C. Each site has a counter, 
initially set to zero, which records the number of arrivals 
on it. We allow the walkers to perform M steps each. 
Then we record the values of each site’s counter which 
are our “daily” fluxes, we reset the counters to zero and 
start again for D “days”. We calculate the average flux 
and standard deviation for each node. Figure shows 
Monte Carlo Simulation results on a scale-free network 
with 7 = 3.0, kmin = 2 and size N = 10000. We allowed 
W = 5000 walkers to simultaneously perform 100 inde¬ 
pendent M = 1000 step walks. We plot the ‘time series’ 
of the number of the walker arrivals at the end of each M- 
step walk (i.e. 100 “days” in this case) at the counters of 
2 sample nodes with degree k = 46, 30 respectively. The 
figure also indicates the average flux (/) and standard 
deviation a for each of the 2 nodes. The main topic of 
this paper consists in the analysis of the fluctuations of 
such ‘time series’. 


Exact results for the star network 

Initially we study the walker dynamics on a simple 
construction such as the star network shown in Fig{^. 
This is a simple case of a bipartite graph with one central 
‘hub’ node (node 0) with degree ko = 9 and 9 other nodes 
with degree k = 1 directly connected to node 0. 

The dynamical problem we have described is exactly 
solvable in the case of the star network. Without loss of 


generality we may assume that the number of steps M is 
an even number. Let as initially consider the case where 
AW = 0 and W = 1, i.e. a single walker performing 
M steps on the star network. In such a case the flux 
of node 0 becomes deterministic because the node will 
be visited in every second step. Thus, if the walker is 
initially placed on any of the nodes 1 to 9, node 0 will be 
visited M/2 times (at steps 2,4,..., M/2). If the walker 
is initially placed on node 0 then the node will again 
be visited M/2 times (at steps 1,3, ...,M/2 — 1). Thus, 
(/o) = M/2 and do = 0 where the averages are over 
different realizations of the M step walk. The fact that 
the dynamics of the central node become deterministic 
(cTo = 0) allow us to calculate the probability that a 
peripheral node is visited m times. The central node is 
accessed exactly M/2 times and each time the peripheral 
node i will be visited with probability p = 1/ko i.e. p = 
1/9 since node 0 has ko = 9. Thus, the probability q{m) 
that a node i, {i = 1,2, ...,9) is visited exactly m times 
is given by a Binomial distribution 

q{m-, M/2,p) = (3) 

The case of W > 1 non-interacting walkers, when the 
network nodes have unlimited capacity C can be viewed 
as equivalent to one walker performing WM steps, lead¬ 
ing to 


q{m;WM/2,p) = -p)WM/2-m (4) 

\ m J 

The mean number of visits on the peripheral nodes and 
the variance are given from the well known formulas 
for the Binomial distribution leading to 


, WM 

fi = ^p 


( 5 ) 


WM 


p{l-p) 


(6) 


Thus, for the case of AW = 0 we recover the power 
law dispersion-flux scaling with exponent a = 1/2 as can 
readily be seen from eqs [5][^ 

For discussing dynamics with “external” noise we 
study the case AW 7^ 0. Let us, for specificity, examine 
the case AW = 1. In this case, the number of walkers 
on our system is a random variable with discrete range 
W — 1, W, W + 1 and probability 1/3 for each of these 
values. 

The resulting distribution for the flux on node i is a 
mixture distribution i.e. the probability distribution of a 
random variable whose values can be interpreted as being 
derived from an underlying set of other random variables 
each with ‘weight’ re = 1/3 in this case. 
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In case of a mixture of one-dimensional normal distri¬ 
butions with weights Wi , means and variances , the 
total mean and variance will be: 

n 

E[X] = fi = '^Wiiii (7) 

i=l 

n 

E[{X - Ilf] = 0-2 = '^Wi{{fii - Ilf + a}) (8) 

i=l 

where E[X] denotes the expected value of random vari¬ 
able X. For the central node we do not have to resolve to 
the above formula because = 0 and we can calculate 
the resulting variance from the variation of walkers in a 
straightforward manner. This calculation yields identical 
results to those obtained with the use of the above eqs|^ 

[6l 

In the following Table [I] we present results for the mean 
flux and dispersion for a set of IF = 5 random walkers 
performing M = 10 step walks on the star network pic¬ 
tured in Fig|^ Exact results were obtained from eqs. [5][8| 
and are in excellent agreement with Monte Carlo simu¬ 
lations of these walks on the star network. 



AIF = 0 

AIF = 1 

Node id 

(/) 


if), 


node 0 

25 

0 

25 

4.08 

nodes 1-9 

2.77 

1.57 

2.77 

1.63 


TABLE I: Mean flux and dispersion for a set of IF = 5 
random walkers performing M = 10 step walks on the 
star network pictured in Fig{^ Numbers are exact 
results using eqs. [Sj^and are in excellent agreement 
with Monte Carlo simulations of these walks on the star 
network. 

These results give us immediately an indication of the 
effect of the “external” noise on the dynamics. First, we 
observe that the high degree node gets the majority of 
the flux. Moreover, we notice that the mean flux does not 
change with the introduction of a non zero AIF. Finally, 
we observe that the impact of AIF on the variance is 
strongly dependent on the degree of the node. The low- 
degree nodes have a very mild increase of their variance 
(from 1.57 to 1.63 in our case-study) while the variance 
of the ‘hub’ node jumps from 0 to 4.08. 

The exact results on the star network give us a hint 
on the origin of the a = 1/2 exponent, which is to be 
expected in processes where the observable random vari¬ 
able follows a Binomial (or Poisson) distribution as seen 
in our case study. We can also understand the origin of 
the a = 1 exponent from the following argument. Let F 
denote the total flux on the network i.e. the sum of the 
flux on all nodes. When AIF = 0, F has a fixed value 


equal to the product IFM, i.e. the total number of steps 
from all walkers on the network. When AIF 7^ 0, F be¬ 
comes a random variable. The total flux is distributed 
over the N network nodes (indexed from 0 to A — 1) 

AT-l 

(^)=E(/*) (9) 

i=0 

We set Bi = ^ Bi < (F). We write Bi = 

ai • (F) with 0 < < 1. Then, 

{fi)=ai{F) (10) 

For example, in our case-study for the ‘hub’ of the star 
network = 1/2 since {fi) = 25 and (F) = 50. This 
holds for every AIF. For the rest of the nodes Eqj^ and 
Table [I| lead to, • 50 = 2.77 ^ = 0.055. When 

AIF = 0, this reduces to {fi) = a^F = aiMW. 

We know, however, that if F, X are random variables 
and Y = aX where a is a constant then the variances of 
the two are connected by ((F^)) = a^((A^)), where we 
use the double bracket notation to denote the variance i.e. 
(jy = = {Y‘^) — Thus, when F is variated 

by introducing AIF > 0 we expect from Eq[l^ 



(11) 

From eqs 10||ll|we obtain: 


(if!)) = 


(12) 

to 

II 

to 

to 

to 

(13) 

^ 


(14) 


leading to the desired scaling form. Eqj^can be read- 
ily verified for the star network using the results of Table 

m 


MONTE CARLO SIMULATIONS AND RESULTS 

The above results give us an indication of the way the 
two commonly observed scaling exponents appear and a 
hint for the mechanism of the transition from the one lim¬ 
iting case to the other. They are a helpful guide for the 
computational study of larger networks. We have seen 
that a = 1/2 scaling is to be expected when our observ¬ 
able quantity F follows a Binomial or Poisson distribu¬ 
tion while a = 1 scaling arises in the case of random vari¬ 
ables F, X having a multiplicative relation F = hX with 
constant a since then (F) = b{X) and ((F^)) = lP‘{{X‘^)) 
which eliminating b will lead to ((F^)) = 
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FIG. 2: Dispersion and flux versus degree.(A) Mean flux 
(/) as a function of the node degree k for networks with 
N = 10^, 7 = 2.0, 3.0 and AVF = 0.(B)Flux standard 
deviation a as a function of the node degree k for 
networks with 7 = 3.0 and AkF = 0,1000, 2000 (C) Flux 
standard deviation a as a function of the node degree k 
for networks with 7 = 2.5 and AkF = 0, 2000.(D) Flux 
standard deviation <j as a function of the node degree k 
for networks with 7 = 2.0 and AkF = 0,100, 2000. 


We have also seen that large-degree nodes are more 
sensitive when AW > 0 i.e. they are more sensitive to ex¬ 
ternal noise. Hence, a plausible assumption on the influ¬ 
ence of network structure on the dynamics is the follow¬ 
ing. Due to the network heterogeneity a broad flux dis¬ 
tribution becomes observable since well-connected nodes 
get more flux than low-connected ones. The presence of 
external noise influences the nodes in a non-uniform way, 
with highly-connected node fluctuations enhanced much 
more, and thus leading to a change of the dispersion-flux 
relation. 

To verify this assumption for larger networks we have 
simulated diffusion on the largest connected component 
of scale-free networks with nodes A = 10^ and 7 = 
2.0, 2.5,3.0. A number of walkers w G [IF — A IF, IF + 
AIF], with IF always chosen equal to half the largest 
connected component size, is initially randomly placed 
on the network. Walkers perform 100 walks of M = 1000 
steps each. We monitor the number of visits on each node 
at the end of the M steps obtaining a series of fluxes for 
each node. We calculate the mean flux (/) and the flux 
standard deviation a for each node. 

FigJ^ shows the mean flux (/) as a function of the 
node degree k. As intuitively expected (/) is proportional 
to k reflecting the fact that nodes with higher connectiv¬ 
ity are more frequently visited. This result is valid for 
all scale-free exponents 7, thus, we observe straight lines 


FIG. 3: (A) Node flux standard deviation a as a 
function of the mean flux (/) for scale free networks 
with N = 10^, 7 = 3.0 and AVF = 0,150,1000. (B) 
Node flux standard deviation a as a function of the 
mean flux (/) for scale free networks with N = 10^, 
7 = 2.0 and AVF = 0,150,1000 


with slope equal to 1 in this double logarithmic plot. It is 
also in agreement with a theoretical derivation presented 
in 117]. Figs(2^ and C show the flux standard deviation 
cr as a function of the node degree k for networks with 
7 = 3.0 and AVF = 0,1000, 2000(top,right) and 7 = 2.5 
and AVF = 0, 2000. We observe that a is well described 
as a power law of the degree k of the node. The observed 
slopes are 1/2 for AVF = 0 and 1 for AVF = 1000,2000. 
This is in accordance with the intuition gained from the 
exact results obtained for the star network and the fact 
that (/) is proportional to k. Fig§D shows the same for 
7 = 2.0 and AVF = 0,100, 2000. In this case, due to the 
broad k distribution one can see the different sensitivity 
of the nodes to external noise. The middle curve shows 
one regime of high-degree nodes that have been consid¬ 
erably affected by AVF and another regime of low-degree 
nodes that remain practically unaffected with a data col¬ 
lapse of the low part of the curves for AVF = 0,100. 
Arrows indicate the considerable ‘shift’ of a of the large 
degree nodes when AVF is increased. 

In Figwe plot the node flux standard deviation a as a 
function of the mean flux (f) for scale free networks with 
N = 10^ and 7 = 3.00 (Figj^), 7 = 2.0 (Figj^) for dif¬ 
ferent values of AVF. We observe an intermediate regime 
of curves with slopes considerably different from 1/2 or 
1 depending on the value of AVF. In the case of 7 = 3.0 
one may consider a single power law with slope 0.66 that 
describes adequately the simulation data of AVF = 150 
for the whole range of (/) with the exception of the very 
low fluxes. For 7 = 2.0 and AVF = 150 there is obvi- 
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FIG. 4: (A) Variance versus flux (/) for networks 
with 7 = 2.0, 2.5, 3.0, kmin = 2 and size N = 10^ 
(squares,circles,triangles), 7 = 3.0, kmin = 2 and size 
V = 10^ (stars) and 7 = 3.0, kmin = d and size N = 10^ 
(cross) for walks with AW = 0. The line is the 
theoretical prediction = (/). (B) Variance versus 
flux (/) for networks with 7 = 2.5, kmin = 2 and size 
N = 10^ for 3 different AW = 10,100,1000. Lines are 
fitting to Eq 0 

ously a cross-over from low (/) data that scale with an 
exponent 1/2 to a regime of high (/) that scale with an 
exponent 0.91. 

Thus, it is obvious that a single power law is not an 
adequate description of the flux-dispersion relation in all 
cases. It may be sufficient, however, when someone is in¬ 
terested in the behavior at the asymptotic limit of large 
fluxes. While an exact result is not available for large 
scale free networks, in contrast to the star network case, 
Eqis a rather plausible alternative as shown in [37] . 
There, the authors approximate the arrivals on a node 
when AW = 0 with a Poisson process. Then the case of 
AW 7 ^ 0 is treated as a mixture distribution and Eqj^or 
the variance as a function of the flux is derived. Eq^ is 
expected to hold when the arrival statistics for AkF = 0 
are not considerably different from that of a Poisson dis¬ 
tribution i.e. for large networks with an adequate num¬ 
ber of steps performed by the walkers. The theoretically 
predicted [37] value of the parameter c is (AW/W)‘^. 

In Eigj^we plot the results of large scale Monte Carlo 
simulations of random walks on scale-free graphs in order 
to examine the range of validity of Eqj^ EigJ^ shows 
versus flux (/) for networks with 7 = 2.0, 2.5,3.0, kmin = 
2 and size N = 10^ (squares,circles,triangles), 7 = 
3.0, kmin = 2 and size N = 10^ (stars) and 7 = 
3.0, kmin = b and size V = 10^ (cross) for walks with 
AW = 0. The line is the theoretical prediction = (/) 


i.e. EqJ^for AW = 0 (c = 0) according to [37]. This 
equality is a direct consequence of the well known prop¬ 
erty of the Poisson distribution to have a variance equal 
to its mean value. 

Note that the plot scale is doubly logarithmic. Al¬ 
though the vertical distance of the points from the line 
may seem small on visual inspection it is actually quite 
significant due to the log scale of the y-axis. 

We see that for networks with 7 = 2.0, kmin = 2 and 
7 = 3.0, kmin = b the points fall on the straight line indi¬ 
cating that = (/) is valid. Eor 7 = 2.5, 3.0, kmin = 2, 
however, the vertical distance from the line is quite differ¬ 
ent from zero even for very large networks with V = 10 ^ 
nodes and one may observe 3(/). The reason for 

this difference is rooted in the discrete nature of the ar¬ 
rival statistics which as can be seen from the star network 
is actually described by a Binomial distribution. The Bi¬ 
nomial distribution is known to coincide to a Poisson 
distribution when the probability of success (in this case 
arrival) tends to zero. Eor networks with large average 
degree (i.e. 7 = 2.0, kmin = 2 and 7 = 3.0, kmin = b) 
the mean arrival probability is small and we see a good 
agreement with the theoretical derivation which is based 
on the assumption of a Poisson distribution of the ar¬ 
rivals. In Eig. 1^ we plot versus flux (/) for networks 
with 7 = 2.5, kmin = 2 and size A = 10^ for 3 different 
AW = 10,100,1000. Lines are fitting to Eqj^ with one 
adjustable parameter, the parameter c. We see that Eqj^ 
is a very effective way of analyzing the data compared to 
a power law fitting (Eq. [^. In any case, we believe that 
using a power law in order to describe such data, as is 
routinely done in the literature is legitimate providing 
one keeps in mind that it is an approximate law and is 
usually a sufficient description only if one is interested 
mainly for the behavior when (/) is large. 

Using such a power law may indeed elucidate some 
characteristic properties of the system. Figj5] shows 
Monte Carlo simulation results for the exponent a of 
the dispersion-flux relation (EqQ as a function of the 
“external” noise parameter AW for networks with 7 = 
2.0, 2.5, 3.0 and minimum degree kmin = 1 (FigJ^) or 
kmin = 2 (Eigj^). The choice kmin = 2 makes the 
largest cluster of the networks equal to the network size 
N. Thus, the mean number of walkers is identical no 
matter the 7 exponent. Also in the case of kmin = 1 we 
have used different network sizes so that their largest con¬ 
nected components Nic have roughly equal sizes. In the 
simulated cases presented in Eigj^ we had Nic ^ 9700 
for 7 = 2.0, Nic - 9500 for 7 = 2.5 and Nic - 10300 
for 7 = 3.0. We observe that the exponent a approaches 
1 with different rates for different 7 exponents and that 
the more heterogeneous networks, i.e. those with lower 7 
exponents amplify the effect of the noise parameter AW. 
This is in accordance to our expectations from the solu¬ 
tion on the star network, since scale-free networks with 
7 = 2.0 may be thought of as (connected) collections of 
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FIG. 5: (A) Exponent a of the dispersion-flux relation 
as a function of AW for networks with 7 = 2 . 0 , 2.5, 3.0 
and minimum degree kmin = 1- (B) The same but for 
networks with minimum degree kmin = 2 

many ‘stars’ while for 7 = 3.0 those stars are fewer and 
smaller. 

The parametrization Eqj^ is rather useful for visual¬ 
ization purposes as it allows for a parameter Vi = 
to be associated with each node. In Eigj^ we plot a 
scale-free network with 1000 nodes 7 = 2 , kmin = 1 and 
AW = 0, 50,100, 500. The color of the nodes has been 
chosen according to the value of each node i. 

Note that in contrast to the exponent a of Eqj^ which 
represents a collective graph property, the quantity Vi sig¬ 
nifies a node property, hence the subscript i. Actually, 
Vi can be interpreted as the ‘effective’ scaling of node i. 
Nodes with Vi < 0.6 are shown as blue, 0.6 < < 0.75 

are green,0.75 < < 0.85 are yellow and > 0.85 are 

red. A ‘spring layout’ algorithm has been used to posi¬ 
tion the network nodes and thus nodes with high degree 
are placed towards the center of the figure. Eor AIT = 0 
is rather homogeneous with all nodes towards the ‘blue’ 
spectrum. Eor AIT = 50 we observe several well con¬ 
nected nodes towards the center to appear with yellow 
color indicating higher values. These are the nodes 
that have been affected by the introduction of the ‘exter¬ 
nal’ noise. The effect becomes much more dramatic for 
AIT = 100 with most of the central nodes having higher 
Vi values and a characteristic pattern of yellow central 
nodes appearing. Eor AIT = 500 all nodes are actually 
affected from external noise, the figure is rather homoge¬ 
neous again but all nodes are towards the ‘red’ spectrum 
values Vi 1.0. 

Since, moreover, a power law is often used to analyze 
flux-dispersion data in the literature it would be instruc¬ 


tive, and potentially useful in practice, to see how the 
parameter c of Eq|^ and the exponent a of Eq{^ are re¬ 
lated. We can derive an equation in closed form as fol¬ 
lows: Let yi = In (/ + c/^) and ^2 = In Then we 

want to estimate the values of b and, most importantly, 
a that minimize the integral of the squared differences 
I = J^(y 2 — yi)‘^df where m is the maximum flux. The 
logarithms of ^1,^2 are taken for convenience, otherwise 
a closed form relation cannot be obtained. The min¬ 
imization condition is, of course, equivalent to setting 
dl/db = 0 and dljdn = 0. 

When IT walkers are performing random walks on a 
network of N nodes , at the stationary regime the average 
number of walkers on a node with degree k is kW/{Nk) 
where k is the mean degree of the nodes. We can use this 
to estimate the maximum flux m as m = kmaxsW/{Nk) 
where s is the number of steps that the walkers perform. 
We find that the exponent a = n /2 of Eqg (the factor 
of 2 comes from the fact that EqU] relates the standard 
deviation to the flux while EqJ^the variance to the flux) 
is related to the parameter c as : 

a =-^(6Li2 ( —i—) + 12cm+ 

12cm \cm + 1 / 

31 n(cm + l)(ln(cm + 1) — 2(ln(c) + ln(m))) — tt^) 

(15) 

which along with m = kmax^W/{Nk) is the desired 
result connecting the 2 different model parameters to 
the network topology. In Eq{^ Li2 is the polylogarith- 
mic function defined by ^ for n = 2 . 

Eigure is a plot of the exponent a versus parameter 
c for networks with 7 = 2 . 0 , 2 . 5 , 3 . 0 ,= 2 and size 
N = 10 ^ (squares,circles,triangles). Points are estimates 
of the parameters from fitting of Eq{^ (exponent a) and 
Eq. (parameter c) to Monte Carlo Simulation data. 
Lines are the predictions of Eqpls) We observe a rather 
good agreement between the two showing that, despite 
it’s rather complicated functional form, Eq. [^is a useful 
tool in connecting the two ways of data analysis. 

Equation \TE\ can be used in order to estimate the effect 
of network size N on the value of the exponent a. The 
parameter c does not depend on the network size. Its 
theoretical value equals the ratio (AIT/IT)^ and is, thus, 
independent of N and our simulation results confirm this 
independence within statistical errors. The exponent a 
is, however, dependent on N since it is a function of the 
maximum node degree kmax through the dependence of 
m to kmax 5 (see Eqj^. Note that, this effect may be 
difficult to observe in moderate network sizes since, for 
example, a 5 -fold increase of N for a network of 7 = 3.0 
will only lead to a 2-fold (or less) increase of kmax- 

Figjs] shows the exponent a as a function of AIT/IT 
for networks of 7 = 3 . 0 ^ kmin = 2 (EigJ^), 7 = 
2 . 5 , kmin — 2 (Figj^) and 7 = 2 . 0 , kmin = 1 
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FIG. 6 : Scale-free network with 1000 nodes 7 = 2 , kmin = 1 for AVF = 0 , 50,100, 500. Edges are depict as red lines 
connecting nodes. The color of the nodes depends on the value Vi = of each node. Nodes with < 0.6 are blue, 

0.6 < Vi < 0.75 are green,0.75 < r < 0.85 are yellow and Vi > 0.85 are red. 
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FIG. 7: Exponent a versus parameter c for networks 
with 7 = 2.0, 2.5, 3.0, k^in = 2 and size N = 10^ 
(squares,circles,triangles). Lines are predictions of Eqjl^ 


(FigJ^) for different network sizes, namely N = 
1000,5000,10000,50000,100000. We observe that for 
7 = 3.0 it is not possible to observe an increased expo¬ 
nent with network size in the intermediate regime. For 


FIG. 8 : (A) Exponent ct as a function of AW/W for 
networks of 7 = 3.0, kmin = 2 for different network sizes, 
namely N = 1000, 5000,10000, 50000,100000. (B)The 
same for networks with 7 = 2.5, kmin = 2 (G) The same 
for networks with 7 = 2 . 0 , kmin = 1 


7 = 2.5 the data collapse is not so strong and such an in¬ 
crease is observed but still statistical errors do not allow 
the extraction of clear results just from the simulation 
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FIG. 9: Flux standard deviation a versus (/) for a 
network with 7 = 3.0 with N = 10 ^, M = 1000. Largest 
cluster {Nic — 7000) is used. Walkers are initially 
placed randomly on the network as long as the capacity 
of a node permits it. 


data. Eq{l5] is a valuable tool in this case showing that 
the expected increases of the exponent are of the order 
of 0.05 which is well within the statistical error for a. 
For 7 = 2.0 one may observe an increase of the exponent 
in the intermediate regime in agreement with the effect 
predicted by Eqim as these networks have rather large 
kmax values even for medium sized networks. This result 
confirms that the intermediate exponent region will not 
vanish for medium size networks and, thus, signifies an 
important regime, potentially observable in the analysis 
of real systems. 

The parameterization assumes a homogeneous 

network, where all nodes obey the same scaling relation¬ 
ship. The parameterization Eqj^ is capable of account¬ 
ing for two groups of nodes with different behaviors. As 
the average flux (/) is proportional to the degree k the 
two regimes visible in Figsl^ show that fluctuations at 
high-degree nodes and at low-degree nodes scale differ¬ 
ently. Can we have more groups of nodes with distinct 
behaviors under random walks in the network? Surpris¬ 
ingly, this is true if we consider the concept of node ca¬ 
pacity. In all the above cases we have assumed that the 
nodes can support an unlimited number of walkers with¬ 
out performance degradation. In most real-life large scale 
transport systems, however, the nodes have a limited ca¬ 
pacity. Computer routers for example have upper limits 
on the amount of incoming and outgoing traffic rate that 
they will support. Thus, we are interested in studying 
the influence of the capacity C on the flux-dissipation 
relation. For this study we again use scale free network 
topologies and we additionally assume that each node 


has a capacity C defined as the maximum number of 
random walkers that can occupy the node at the same 
time. When C = 1 we have the well-known case of ex¬ 
cluded volume interactions. For the simulations we have 
used the largest cluster of a 10000-node scale free net¬ 
work with 7 = 3.0. The size of the largest cluster is 
Nic ^ 7000. The number of walkers placed on the net¬ 
work was set equal to Nic and AkF = 0. The maximum 
degree kmax of fho network is kmax — 150. Walkers (try 
to) perform M = 1000 steps each. At the end of the steps 
fluxes are recorded for each node. The process repeats 
for 100 times and time-average fluxes and standard devi¬ 
ations are calculated for each node. In Figwe plot the 
flux standard deviation a versus (/) for a network with 
7 = 3.0 with N = 10^, M = 1000. Walkers are initially 
placed randomly on the network as long as the capacity 
of a node permits it. 

We find 3 interesting regimes: 

(a) For C ^ iW/N)kmin the observed flux range is 
decreased but the power law scaling with exponent 1/2 
remains. In this case, (see C = 2 in FigJ^ all nodes feel 
the limitations of the capacity C. Walkers will often try 
to move but will find the destination to be fully occupied. 
In such a case they will not perform a step, resulting in 
a reduced number of arrivals on all counters. 

(b) For C ^ {W/N)kmax we notice the appearance 
of outliers at the beginning and end of the series (see 
C = 50,100 in FigJ^ which do not allow to claim that 
the scaling relation is well approximated by a power-law. 
Here the high degree nodes (hubs) are influenced from the 
capacity limitation. Thus, while the low-degree nodes re¬ 
ceive more or less the same amount of walker arrivals as 
in the unlimited case the hubs reach the capacity limit. 
This limitation significantly alters the visitation pattern 
and leads to the observed appearance of the outliers. In 
several cases the walkers will be unable to visit the sat¬ 
urated hubs and thus will remain unmovable for some 
steps leading to decreased “flux” recorded at the coun¬ 
ters compared to the unsaturated case. 

(c) For C ^ {W/N)kmax Ihe power law scaling with 
exponent 1/2 is, as expected, recovered since the nodes 
can accommodate all possible arrivals without problem 
similar to the case of unlimited capacity. 

Thus, when the capacity parameter is taken into ac¬ 
count, we can roughly distinguish 3 types of nodes i.e. 
high degree saturated nodes, high-medium degree unsat¬ 
urated nodes and low degree nodes, which influence the 
fluctuations scaling in quantitatively different ways. In 
Figure [Tq| we plot 2 extreme capacity cases of a scale-free 
network with 1000 nodes 7 = 2 , kmin = 1 for AW = 0 
and capacity C = 1000(Figs [l^,C), C = 10(Figs 
[TQ^,D)) . High degree nodes are placed near the center of 
the graph. The color of the nodes depends on the value 
Vi = of each node i. Nodes with < 0.45 are yellow, 
0.45 < Vi < 0.65 are blue,and those with > 0.65 are 
red. The case C = 1000 (FigpT|4) is actually a network 
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C=1000 C=10 



FIG. 10 : Scale-free network with 1000 nodes 
7 = 2 , kmin = 1 for AVF = 0 and capacity C = 1000 (A), 
C = 10 (B) . Edges are depict as red lines connecting 
nodes. The color of the nodes depends on the value 
Vi = of each node. Nodes with < 0.45 are yellow, 
0.45 < Vi < 0.65 are blue,and those with > 0.65 are 

red.(C) Histogram of the values for C = 1000 .(D) 
Histogram of the values for C = 10. 


with internal fluctuations only, since AIF = 0, and prac¬ 
tically unlimited capacity. As expected, the fraction is 
close to 0.5 for all nodes and all nodes appear with blue 
color. When C = 10 (Fig[TQ^), although AkF = 0 and no 
external fluctuations are present, we observe the appear¬ 
ance of highly connected saturated nodes (yellow) with 
Vi < 0.45 coexisting with unsaturated nodes (blue nodes 
towards the center of the graph) as well as some nodes 
with high Ti ratio consistent with our previous remarks. 
The bottom figures are histograms of the values for 
the two cases C = 1000(Fig[^) and C = 10 (Fig[l^). 
Notice the appearance of additional ‘bands’ close to 0.25 
and 0.75 which are not present in the case of ‘unlimited’ 
capacity. 

In case we include some external noise AkF in addi¬ 
tion to the capacity we expect to observe a combina¬ 
tion of slope changing due to external noise and the ap¬ 
pearance of ‘outliers’ i.e. unconventional points of the 
cr — f plot due to the saturated nodes. In Figure pj] we 
plot the (7 versus (/) for a network with 7 = 3.0 with 
N = 10 ^, M = 1000 for different capacity C and AkF 
combinations, namely C = 100 , AVF = 100 (small dia¬ 
monds), C = 100 , AkF = 500(large diamonds), C = 
50,A1F = lOO(circles), C = 100, AkF = 500(stars). 
We can again observe the influence of the external noise 
which is increasing the slope of the curves with increasing 
AkF. Especially for the cases of AW = 500 we observe 
that for (/) < 10^ there is a regime with slope around 0.7 



FIG. 11: Flux standard deviation a versus (/) for a 
network with 7 = 3.0 with N = 10^, M = 1000. 
Different capacity C and AW combinations are shown. 

(dashed line) due to the influence of the external noise 
on the low degree nodes, a second regime with a higher 
slope due to the influence of the external noise on the 
high degree unsaturated nodes and out-liner points due 
to the high degree saturated nodes. 

CONCLUSIONS 

Stochastic processes on networks exhibit fluctuations 
due to combinations of internal and external noise. We 
have used a multiple random walk model to study the 
effect of network heterogeneity on the fluctuations of net¬ 
work dynamics and used random walks as a tool for un¬ 
derstanding the relationship between topology and dy¬ 
namics. We have obtained exact results for the star 
network which are indicative of the behavior of large 
scale-free networks. We have found that the network 
heterogeneity acts as an amplifier of the effects of ex¬ 
ternal noise. These effects include a range of exponents 
between 1/2 and 1 and are persistent for medium size 
networks. Moreover, we propose that an analysis of the 
flux variance as a sum of 2 power laws .i.e. a relation 
of the form = (/) + c((/))^ with only one adjustable 
parameter c can provide a more adequate description of 
the problem under investigation. In particular, we de¬ 
rive a semi-analytical expression for the error of ‘data’ 
following equation when describing them with EqUl 
In this way we can understand (i) why the ‘transition’ 
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from Qf = l/ 2 toQf = lat increasing AW is not becom¬ 
ing sharper with increasing network size, and (ii) why 
network heterogeneity amplifies the influence of external 
noise. Finally, we have shown that when the internal 
dynamics correlate with the external influence, as in the 
case of nodes with a maximum capacity, there appear 
regimes with non-power law scaling characterized by the 
appearance of outliers. 

An alternative implementation of a finite capacity is 
to allow queues to form at nodes in times, where the 
total capacity of the node has been exceeded. Such an 
extension of random walks to queuing has been discussed 
in [31]. 

In that case, one can also resort to various analytical 
results in queueing theory. For example, an interesting 
alternative treatment of the scaling crossover in the case 
in which network dynamics is best modeled in terms of 
a network of waiting lines can be found in [JH] . There, 
a system of M departments is considered, each one with 
it’s own mean arrival rate and mean serving time. Cus¬ 
tomers from outside the system arrive at a department 
k in di poisson type manner and when served are moved 
(instantaneously) to another department m with a cer¬ 
tain probability 0(k^m) or the service is completed. The 
steady-state queue length of such a system is shown to 
follow a negative binomial (Pascal) distribution. Assum¬ 
ing that the flux is proportional to the queue length, the 
flux-fluctuations relation can be controlled by changing 
0{k,m) from 0 (non-interacting departments) to 1 (fully 
interacting departments). For non-interacting depart¬ 
ments the arrivals are a poisson process, hence leading 
to exponent 1/2 while the fully interacting case is ex¬ 
pected to lead to exponent 1 for appropriate choices of 
the arrival rates and serving times. 

Exploring scaling relations at the intersection of these 
two topics, queing theory and random walks on graphs, 
further can lead to additional interesting insights into 
the ‘patterns’ on a graph formed by the different scaling 
behaviors. 


[ 1 ] R. Albert and A.-L. Barabasi, Reviews of modern physics 
74, 47 (2002). 

[ 2 ] A. Arenas, A. Dfaz-Guilera, J. Kurths, Y. Moreno, and 

C. Zhou, Physics Reports 469, 93 (2008). 

[3] A.-L. Barabasi, Linked: The New Seienee Of Networks 
Seienee Of Networks (Basic Books, 2002 ). 

[4] S. Boccaletti, V. Latora, Y. Moreno, M. Chavez, and D.- 
U. Hwang, Physics reports 424, 175 (2006). 

[5] S. Bornholdt, H. G. Schuster, and J. Wiley, Handbook of 
graphs and networks, vol. 2 (Wiley Online Library, 2003). 

[ 6 ] M. Newman, Networks: an introduetion (Oxford Univer¬ 
sity Press, 2010). 

[71 M. Bogachev and A. Bunde, Europhysics Letters 86, 
66002 (2009). 

[ 8 ] L. Chen, J. Chen, Z.-H. Guan, X.-H. Zhang, and 


D. -X. Zhang, Physica A: Statistical Mechanics and 

its Applications 391, 3336 (2012), ISSN 0378- 

4371, URL http://www.sciencedirect.com/science/ 
article/pii/S0378437112000052 

[9] B. Tadic, G. Rodgers, and S. Thurner, International 
Journal of Bifurcation and Chaos 17, 2363 (2007). 

[10] D. Helbing, D. Armbruster, A. S. Mikhailov, and 

E. Lefeber, Physica A: Statistical Mechanics and 
its Applications 363, xi (2006), ISSN 0378-4371, 
information and Material Flows in Complex Net¬ 
works, URL http://www.sciencedirect.com/science/ 
article/pii/S0378437106000835 

[11] D. Helbing, A. Deutsch, S. Diez, K. Peters, 
Y. Kalaidzidis, K. Padberg-Gehle, S. Lammer, A. Jo¬ 
hansson, G. Breier, F. Schulze, et al.. Advances in 
Complex Systems 12, 533 (2009), ISSN 02195259, URL 
http://search.ebscohost.com/login.aspx?direct= 
true&db=aph&AN=47167240&site=ehost-live 

[12] F. Matthaus, C. Salazar, and O. Ebenhoh, PLoS Com¬ 
putational Biology 4, el000049 (2008). 

[13] N. Sonnenschein, C. Marr, and M.-T. Hiitt, Metabolites 
2(3), 632 (2012). 

[14] R. D. Smith, Advances in Complex Systems (ACS) 14, 
905 (2011). 

[15] A. Vazquez, R. Pastor-Satorras, and A. Vespignani, 
Physical Review E 65, 066130 (2002). 

[16] T. Becker, M. E. Beber, K. Windt, M.-T. Hiitt, and 
D. Helbing, Journal of Statistical Mechanics: Theory and 
Experiment 2011, P05004 (2011), URL http://stacks. 
iop.org/1742-5468/2011/i=05/a=P05004 

[17] P. Holme, Advances in Complex Systems 6, 163 (2003). 

[18] A. Vespignani, Science 325, 425 (2009). 

[19] C. Fretter, L. Krumov, K. Weihe, M. Miiller-Hannemann, 
and M. Hiitt, European Physical Journal B: Condensed 
Matter Physics 77, 281 (2010). 

[20] T. Becker, M. E. Beber, K. Windt, and M.-T. Hiitt, 

CIRP Journal of Manufacturing Science and Technol¬ 
ogy 5, 309 (2012), ISSN 1755-5817, special issue 

from 44th CIRP Conference on Manufacturing Sys¬ 
tems, URL http://www.sciencedirect.com/science/ 
article/pii/S1755581712000570 

[21] E. Almaas, B. Kovacs, T. Vicsek, Z. N. Oltvai, and A.-L. 
Barabasi, Nature 427, 839 (2004). 

[22] M. Rosvall and C. Bergstrom, Proc. Natl. Acad. Sci. U. 
S. A. 105, 1118 (2008). 

[23] L. Schulman and B. Gaveau, Bull. Sci. Math. 129, 631 
(2005). 

[24] M. Miiller-Linow, C. C. Hilgetag, and M.-T. Hiitt, PLoS 
Computational Biology 4, el000190 (2008), URL http: 
//dx. doi . org/10.137l7o2Fj ournal. pcbi . 1000190 

[25] D. Brockmann and D. Helbing, Science 342, 1337 (2013). 

[26] J. S. Edwards, R. Ibarra, and B. 0. Palsson, Nature 
Biotechnology 19, 125 (2001). 

[27] N. D. Price, J. L. Reed, and B. 0. Palsson, Nature Re¬ 
views. Microbiology 2 , 886 (2004). 

[28] N. C. Duarte, S. A. Becker, N. Jamshidi, 1. Thiele, M. L. 
Mo, T. D. Vo, R. Srivas, and B. 0. Palsson, Proc. Natl. 
Acad. Sci. U. S. A. 104, 1777 (2007). 

[29] P. Kaluza, M. Vingron, and A. Mikhailov, Chaos 18, 
026113 (2008). 

[30] M. Beber, D. Armbruster, and M.-T. Hiitt, European 
Physical Journal B 86, 473 (2013). 

[31] J. Duch and A. Arenas, Physical Review Letters 96, 
218702 (2006). 



12 


[32] L Simonsen, K. Astmp Eriksen, S. Maslov, and K. Snep- 
pen, Physica A: Statistical Mechanics and its Applica¬ 
tions 336, 163 (2004). 

[33] L. Lovasz, Combinatorics, Paul erdos is eighty 2, 1 
(1993). 

[34] M. A. De Menezes and A.-L. Barabasi, Physical Review 
Letters 92, 028701 (2004). 

[35] Z. Eisler, I. Bartos, and J. Kertesz, Advances in Physics 
57, 89 (2008). 

[36] L. Prignano, Y. Moreno, and A. Diaz-Guilera, Physical 
Review E 86 , 066116 (2012). 

[37] S. Meloni, J. Gomez-Gardenes, V. Latora, and 
Y. Moreno, Physical Review Letters 100, 208701 (2008). 

[38] M. P. Stump! and M. A. Porter, Science 335, 665 (2012). 

[39] E. Almaas and A.-L. Barabasi, in Power laws, scale-free 
networks and genome biology (Springer, 2006), pp. 1-11. 

[40] A. Arenas, L. Danon, A. Diaz-Guilera, P. M. Gleiser, 


and R. Guimera, The European Physical Journal B- 
Condensed Matter and Complex Systems 38, 373 (2004). 

[41] M. Barthelemy, J.-P. Nadal, and H. Berestycki, Proc. 
Natl. Acad. Sci. U. S. A. 107, 7629 (2010). 

[42] M.-B. Hu, W.-X. Wang, R. Jiang, Q.-S. Wu, and Y.-H. 
Wu, EPL (Europhysics Letters) 79, 14003 (2007). 

[43] M.-T. Hiitt and U. Liittge, Physica A: Statistical Me¬ 
chanics and its Applications 350, 207 (2005). 

[44] M. Kampf, S. Tismer, J. W. Kantelhardt, and L. Much- 
nik, Physica A: Statistical Mechanics and its Applica¬ 
tions 391, 6101 (2012). 

[45] V. Kishore, M. Santhanam, and R. Amritkar, Physical 
Review Letters 106, 188701 (2011). 

[46] Y. Xia and D. Hill, Europhysics Letters 89, 58004 (2010). 

[47] J. D. Noh and H. Rieger, Physical Review Letters 92, 
118701 (2004). 

[48] J. R. Jackson, Operations research 5, 518 (1957). 



